Genome-Wide Association Study Identiﬁes Genetic Polymorphisms Associated with Estimated Minimum Effective Concentration of Fentanyl in Patients Undergoing Laparoscopic-Assisted Colectomy

: Sensitivity to opioids varies widely among individuals. To identify potential candidate single-nucleotide polymorphisms (SNPs) that may signiﬁcantly contribute to individual differences in the minimum effective concentration (MEC) of an opioid, fentanyl, we conducted a three-stage genome-wide association study (GWAS) using whole-genome genotyping arrays in 350 patients who underwent laparoscopic-assisted colectomy. To estimate the MEC of fentanyl, plasma and effect-site concentrations of fentanyl over the 24 h postoperative period were estimated with a pharmacokinetic simulation model based on initial bolus doses and subsequent patient-controlled analgesia doses of fentanyl. Plasma and effect-site MECs of fentanyl were indicated by fentanyl concentrations, estimated immediately before each patient-controlled analgesia dose. The GWAS revealed that an intergenic SNP, rs966775, that mapped to 5p13 had signiﬁcant associations with the plasma MEC averaged over the 6 h postoperative period and the effect-site MEC averaged over the 12 h postoperative period. The minor G allele of rs966775 was associated with increases in these MECs of fentanyl. The nearest protein-coding gene around this SNP was DRD1 , encoding the dopamine D 1 receptor. In the gene-based analysis, the association was signiﬁcant for the SERP2 gene in the dominant model. Our ﬁndings provide valuable information for personalized pain treatment after laparoscopic-assisted colectomy.


Introduction
Opioids, such as morphine, fentanyl, oxycodone, and hydromorphone, are widely used as effective analgesics for the treatment of acute and chronic pain because of their robust antinociceptive effects. However, effects of opioids are not uniform across all patients, and considerable differences in the responsiveness or sensitivity to opioids are widely known [1,2]. This can influence the total amount of analgesics that are required for adequate pain relief, which can hamper the effective treatment of pain in clinical practice. For example, the minimum effective concentrations (MECs) of fentanyl at which patients demand additional fentanyl doses to relieve recurring postoperative pain are reported to vary widely among patients from 0.23 to 0.99 ng/mL after orthopedic surgery [3], from 0.23 to 1.18 ng/mL after open abdominal surgery [4], from 0.30 to 1.45 (5-95 percentiles) ng/mL after open abdominal surgery [5], and from 0.2 to 8.0 ng/mL after various surgical procedures [6], indicating that MECs of fentanyl after certain surgical procedures have a more than four-to fivefold difference among individuals [7]. Because of such significant differences in opioid sensitivity, empirical methods of administration that have been utilized by trial and error are an imperfect practice that can result in delayed or inadequate analgesia and possibly overdose [7].
The required amounts of opioid analgesics may also vary among patients with pain depending on age, sex, weight, basal pain sensitivity, the type of surgery, perceived pain during the perioperative period [2], and genetic factors. Previous twin studies of experimental heat and cold pressor pain reported that genetic effects were estimated to account for 12%, 60%, and 30% of the observed response variance (i.e., pain threshold) after administration of the opioid analgesic alfentanil for heat pain, cold-pressor pain, and cold-pressor pain, respectively [8,9]. Although the variance of responses to opioids appears to be moderately influenced by genetic factors, potential genes and genetic variants that are involved in response variance have not yet been fully elucidated. Further studies are needed to delineate such genetic factors.
To date, many candidate gene association studies have been conducted [10][11][12]. These studies have targeted various genes that are involved in pharmacokinetic and pharmacodynamic opioidergic pathways and pain-related genes of various modalities, such as the µ-opioid receptor (OPRM1) gene; cytochrome P450, family 2, subfamily D, polypeptide 6 (CYP2D6) gene; adenosine triphosphate-binding cassette (ABC), subfamily B (MDR/TAP), member 1 (ABCB1) gene; catechol-O-methyltransferase (COMT) gene; and genes that are related to cytokines (e.g., interleukin-1β, interleukin-6, and tumor necrosis factor-α) [2]. Additional candidate genes are detailed in previous review articles [10][11][12]. Genetic factors that are related to opioid sensitivity and responsiveness can also be explored using a genome-wide approach in genome-wide association studies (GWASs), although only a few studies have conducted GWASs of such phenotypes. One example is a prospective cross-sectional multinational multicenter study of patients with cancer from 11 European countries [13] who were treated with opioids for moderate or severe pain. The strongest association with responsiveness to opioids was found for the rs12948783 single-nucleotide polymorphism (SNP), which is located upstream of the RHBDF2 gene [14].
We also conducted GWASs of phenotypes that are related to opioid sensitivity and candidate gene studies [15][16][17][18][19][20][21][22][23][24][25][26][27][28][29][30][31][32][33]. In our GWASs, although genetic variants that were significantly associated with opioid responsiveness for the treatment were not found in patients with chronic pain [30], we identified several SNPs, including the rs2952768 SNP (located near the METTL21A [FAM119A] and CREB1 gene regions), that were significantly associated with postoperative opioid analgesic requirements in subjects who underwent cosmetic orthognathic surgery for mandibular prognathism [18]. Furthermore, a GWAS of patients who were treated with opioid analgesics for the treatment of cancer pain identified several SNPs that were significantly associated with average daily opioid requirements for the treatment of pain, the best candidates of which were the rs1283671 and rs1283720 SNPs in the ANGPT1 gene region. We also conducted GWASs of subjects who underwent laparoscopic-assisted colectomy (LAC), a surgery that is often categorized as minimally invasive because of much smaller skin incisions and less postoperative pain compared with traditional open abdominal surgery, although postoperative pain is not "minimal" after surgery [22]. Our GWASs of subjects who underwent LAC identified several potent SNPs, including the nonsynonymous rs2076222 SNP in the LAMB3 gene region, the rs199670311 nonsynonymous SNP in the TMEM8A gene region, and intronic SNPs, including rs4839603, in the SLC9A9 gene region [22,25].
Likely because of relative facileness, most previous human genetic studies have focused on opioid analgesic requirements for the treatment of disease-related pain, chronic pain, and perioperative/postoperative pain as the main endpoint to investigate genetic variants that are associated with human responsiveness and sensitivity to opioids . However, MECs at the plasma and effect site are not generally measured directly. Thus, these parameters have not been used to date in human genetic association studies. Nevertheless, with the development of pharmacokinetic/pharmacodynamic knowledge and the advancement of computer technology, it has become easier to simulate the process of plasma or effect-site concentrations of anesthetics and analgesics by leveraging related simulation software programs, such as STANPUMP (http://opentci.org/code/stanpump; accessed on 25 January 2023) and tivatrainer (https://www.tivatrainer.com; accessed on 25 January 2023). Such simulation software has been used in many studies to estimate plasma and effect-site concentrations of anesthetics and analgesics [34][35][36][37][38][39].
In the present study, we conducted a GWAS of patients who underwent LAC to identify potential genetic variants that contribute to the efficacy of opioid analgesics based on information about estimated plasma or effect-site concentrations of fentanyl, which were calculated by utilizing one of the programs, the BeConSim Monitoring simulation software program (http://www.masuinet.com; accessed on 1 January 2020) [40][41][42].

Impact of Clinical Variables on Estimated MEC of Fentanyl in Subjects Who Underwent LAC
All 351 subjects completed the study. However, data were incomplete for one subject, particularly postoperative data. Therefore, postoperative data from the remaining 350 subjects were analyzed. Demographic, anesthetic, and surgical data for all 351 subjects are detailed in Supplementary Table S1 and our previous reports [22,25].

Identification of Genetic Polymorphisms Associated with Estimated MEC of Fentanyl in Patients Who Underwent LAC by GWAS
We then explored the association between genetic variations and opioid sensitivity, which was evaluated as the estimated plasma and effect-site MECs after surgery in a total of 350 subjects who underwent LAC that involved the administration of opioid analgesics [22,25]. The surgical procedure was relatively uniform; thus, invasiveness and the resultant pain were regarded as relatively homogeneous among subjects. GWASs were conducted as a consecutive three-stage analysis to identify potent SNPs that were associated with the estimated 0-6 h plasma MEC, 0-12 h plasma MEC, 0-6 h effect-site MEC, and 0-12 h effect-site MEC. Consequently, 14, 26, and 10 SNPs were selected as the top candidates in the additive, dominant, and recessive models, respectively, after the final stage for the 0-6 h plasma MEC (Supplementary Figure S1A). For the 0-12 h effect-site MEC, 14, 28, and 8 SNPs were selected as the top candidates in the additive, dominant, and recessive models, respectively, after the final stage (Supplementary Figure S1B). Similarly, 21, 24, and 19 SNPs were initially selected in the additive, dominant, and recessive models, respectively, for the 0-12 h plasma MEC. Likewise, 17, 25, and 10 SNPs were initially selected in the additive, dominant, and recessive models, respectively, for the 0-6 h effect-site MEC (details not shown). The potent SNP lists are presented in Tables 1 and 2 and Supplementary Tables S2 and S3. Among these, one SNP, rs966775, that mapped to 5p13 (GRCh37) showed significant associations with the 0-6 h plasma MEC and 0-12 h effect-site MEC after the final stage in the additive model (combined β = 0.0916, nominal p = 1.027 × 10 −7 , for the 0-6 h plasma MEC; combined β = 0.1071, nominal p = 1.299 × 10 −7 , for the 0-12 h effect-site MEC; Tables 1 and 2). The observed p values for this SNP, calculated as −log10 (p value), deviated from the expected values from the null hypothesis of uniform distribution in the quantile-quantile (QQ) p-value plots for the entire sample (Supplementary Figure S1 for the 0-6 h plasma MEC, Supplementary Figure S2 for the 0-12 h effect-site MEC). Similar strong associations with this SNP were observed for the 0-12 h plasma MEC and 0-6 h effect-site MEC after the final stage in the additive model (combined β = 0.0908, nominal p = 1.206 × 10 −7 , for the 0-12 h plasma MEC; combined β = 0.1095, nominal p = 1.942 × 10 −7 , for the 0-6 h effect-site MEC; Supplementary Tables S2 and S3), although the associations were not significant. The rs966775 SNP is located in the intergenic region, and the nearest protein-coding gene from this SNP position was DRD1, which encodes the dopamine D 1 receptor (Supplementary Figure S3). A linkage disequilibrium (LD) block that includes the rs966775 SNP was assumed to span the approximately 1 kbp chromosomal region, and no SNPs showed high LD with this SNP in the neighboring region, including the DRD1 gene region (pairwise calculated r 2 = 0.93; Supplementary Figure S3). When MECs (in ng/mL) were log-transformed and shown as mean ± standard error of the mean (SEM), 0-6 h plasma MECs were 0.4960 ± 0.0171, 0.5289 ± 0.0192, and 0.6839 ± 0.0393, and 0-12 h effect-site MECs were 0.5393 ± 0.0188, 0.5839 ± 0.0220, and 0.7589 ± 0.0433, in subjects with the A/A (n = 171), A/G (n = 137), and G/G (n = 41) genotypes of this SNP, respectively. The copy number of the minor G allele was associated with higher 0-6 h plasma and 0-12 h effect-site MECs. A similar trend was observed for 0-12 h plasma and 0-6 h effect-site MECs, and the copy number of the minor G allele was associated with greater MEC values for these phenotypes. The genotype distribution of this SNP met the criteria of the Hardy-Weinberg equilibrium tests (χ 2 = 2.7283, p = 0.0986).

Identification of Genes and Gene Sets Associated with Estimated MEC of Fentanyl in Patients Who Underwent LAC by Gene-Based and Gene-Set Analyses
Considering that effects of individual markers tend to be too weak to be detected by comprehensive analyses, such as GWASs, which target only single polymorphisms, we conducted gene-based and gene-set analyses, which are statistical methods that are used to analyze multiple genetic markers simultaneously to determine their joint effect. In both analyses using MAGMA software [43], which was made accessible in the FUMA GWAS platform [44], we investigated genes and gene sets that were related to the estimated MEC of fentanyl in a total of 350 patients who underwent LAC. As a result, 921,239 SNPs from the selected candidate genes and gene sets in the additive, dominant, and recessive models were included in the analyses of all patients. The top 20 candidate genes that were found in each genetic model by the gene-based analysis are listed in Table 3. In the dominant model, SERP2, the top candidate gene, was significantly associated with the 0-6 h plasma MEC (adjusted p = 0.02425; Table 3, Figure 1B Table 4, Figure 2B). The association between the SERP2 gene and the 0-6 h effect-site MEC was marginally significant (adjusted p = 0.05245; Supplementary Table S5). However, in both the additive and recessive models, none of the genes were significantly associated with the phenotypes (Tables 3 and 4; Supplementary  Tables S4 and S5; and Figures 1A,C and 2A,C). The top 20 candidate gene sets for each phenotype that were found in each genetic model by the gene-set analysis are listed in Supplementary Tables S6-S9. As a result, the 0-6 h plasma MEC was significantly associated with the "go_paracrine_signaling" (adjusted p = 0.01093), "reac-tome_free_fatty_acid_receptors" (adjusted p = 0.01430), and "go_taste_receptor_activity" (adjusted p = 0.03004) gene sets in the additive model (Supplementary Table S6) and the "go_negative_regulation_of_epidermal_cell_differentiation" (adjusted p = 0.01728), "go_paracrine_signaling" (adjusted p = 0.02210), "go_negative_regulation_of_epidermis_ development" (adjusted p = 0.03244), and "go_negative_regulation_of_keratinocyte_ differentiation" (adjusted p = 0.04440) gene sets in the recessive model (Supplementary Table S6). The 0-12 h plasma MEC was significantly associated with the "sotiriou_breast_ cancer_grade_1_vs_3_dn" gene set in the additive model (adjusted p = 0.04805; Supplementary Table S7). The 0-6 h effect-site MEC was significantly associated with the "go_ccr2_chemokine_receptor_binding" and "go_paracrine_signaling" gene sets in the additive model (adjusted p = 0.01342 and 0.03030, respectively; Supplementary Table S8) and significantly associated with the "go_paracrine_signaling" and "go_ccr2_chemokine_ receptor_binding" gene sets in the recessive model (adjusted p = 0.01030 and 0.02499, respectively; Supplementary Table S8). The 0-12 h effect-site MEC was significantly associated with the "pid_shp2_pathway" gene set in the recessive model (adjusted p = 0.02854; Supplementary Table S9). The genes that were included in these gene sets are listed in Supplementary Table S10. The SERP2 gene, which was significantly associated with the phenotypes in the gene-based analysis, was not included in any of the gene sets (Supplementary Table S10). Among these genes, several genes were commonly included in two or three kinds of gene sets (Supplementary Table S10). Eight genes (EZH2, GRHL2, HOXA7, MSX2, REG3A, REG3G, SRSF6, and TP63) were commonly included in three kinds of gene sets (Supplementary Table S10).

Discussion
Although human genetic variants that are associated with human responsiveness and sensitivity to opioids have been explored by adopting opioid analgesics that are required for the treatment of disease-related pain, chronic pain, and perioperative/postoperative pain as the main endpoint , plasma and effect-site MECs have not been investigated in genetic studies, likely because of difficulties in measuring actual values of plasma and effect-site MECs. To comprehensively explore genetic factors that underlie large individual differences in fentanyl responsiveness and sensitivity after LAC, we first conducted a GWAS in this cohort of patients, focusing on plasma and effect-site MECs of fentanyl that were estimated with a pharmacokinetic simulation model [40][41][42]. As a result of the GWAS in surgical patients, 8-28 SNPs were selected as the top candidate SNPs that were significantly associated with a plasma or effect-site MEC that was averaged over the 0-6 h or 0-12 h postoperative period after LAC in all of the additive, dominant, and recessive models (Tables 1 and 2; Supplementary Tables S2 and S3). Among these, the rs966775 SNP that mapped to 5p13 had highly significant associations with 0-6 h plasma and 0-12 h effect-site MECs (Tables 1 and 2). A gene that is located near the region of this SNP was DRD1, which encodes the dopamine D 1 receptor. Altogether, our data indicated that the rs966775 SNP near the DRD1 gene significantly affected fentanyl sensitivity. Compared with non-carriers, G-allele carriers of this SNP were associated with higher plasma and/or effect-site MECs of fentanyl, suggesting that G allele carriers would feel pain at a higher plasma/effect-site fentanyl concentration and thus would require more frequent self-dosing of fentanyl for adequate pain control. Although we acknowledge that the sample size of 350 patients may not be sufficiently large to draw definitive conclusions about genetic markers that contribute to individual differences in the MEC of fentanyl and that further research is needed with larger sample sizes and greater statistical power to validate our findings, the present results suggest that the rs966775 SNP could serve as a marker that predicts the efficacy of opioid analgesics for the treatment of postoperative pain.
In clinical postoperative pain management using patient-controlled analgesia (PCA), continuous pain relief should be achieved if the plasma opioid concentration is maintained in excess of the MEC, whereas pain will return if it decreases to the MEC. Thus, the MEC is indicated by the need for an additional intravenous (i.v.) opioid because of recurring pain [4,5]. The MEC of opioids varies depending on the type of surgery and intensity of postoperative pain. It gradually decreases with a decreasing intensity of postoperative pain [3][4][5][6]. Nevertheless, the MEC remains relatively constant within each patient over the postoperative period but varies widely among patients even after the same type of surgery [3][4][5]. Associations between genetic variants and MECs of opioids have not been investigated in genetic studies, likely because of difficulties in repeatedly measuring actual plasma opioid concentrations. However, opioids act on the effect site and not on plasma, and pharmacokinetic simulation models can predict plasma concentrations with acceptable accuracy [35,45,46] and estimate effect-site concentrations that are not measurable in humans [34,36,37]. Therefore, simulation models have been widely used in clinical studies [6,[34][35][36][37][38][39], including one that evaluated the plasma MEC of fentanyl [6]. Using MECs of fentanyl that were determined with a simulation model, we conducted a GWAS and found that the rs966775 SNP significantly affected fentanyl sensitivity.
As mentioned above, we conducted GWASs of phenotypes that are related to opioid sensitivity and candidate gene studies [15][16][17][18][19][20][21][22][23][24][25][26][27][28][29][30][31][32][33]. Several candidate SNPs were found to be associated with phenotypes that are related to opioid sensitivity/pain. Among these are the rs2076222 SNP in the LAMB3 gene region, which was associated with postoperative 24 h fentanyl requirements in subjects who underwent LAC. Although this SNP was also found to be nominally significantly associated with the estimated plasma and effect- The associations between other candidate SNPs that we identified in previous studies as candidates for opioid sensitivity and the estimated plasma and effect-site MECs in the present study were not even nominally significant (details not shown). These results might indicate the general difficulty in replicating results of human genetic association studies or reflect phenotypical differences between postoperative analgesic requirements per se and the plasma and effect-site MECs that were estimated with a pharmacokinetic simulation model in the present study.
The best candidate SNP in the present study was rs966775, which was found in the intergenic region. The protein-coding gene on chromosome 5 that was nearest to this SNP site was the DRD1 gene. This gene encodes the dopamine D 1 receptor, which is the most abundant dopamine receptor in the central nervous system. Although the D 1 receptor has been shown to be involved in mechanisms of opioid analgesia in animal studies [47][48][49][50], the impact of this SNP on the expression and function of the DRD1 gene product is not known but presumably may not be profound because this SNP is located more than 100 kbp from the gene region (Supplementary Figure S3). The rs966775 SNP has not been previously reported to be associated with any phenotypes to date. Although this SNP was not in strong LD (r 2 ≥ 0.80) with any other neighboring SNPs in our data (Supplementary Figure S3), when these SNPs were referenced in HaploReg v. 4.1 and SNPinfo Web Server (accessed on 30 January 2023) [51,52], they were in strong LD with the rs7725278, rs897747, rs2382021, rs2890873, rs3955076, rs76895738, and rs10060502 SNPs (r 2 ≥ 0.80) and were moderately linked to the rs12652255 SNP (r 2 = 0.68) in Asian populations, including Japan. HaploReg v. 4.1 also showed that the rs966775 SNP could change six motifs for DNA-binding proteins and overlaps with an enhancer in the fat and skin. Nevertheless, none of these SNPs were significantly associated with mRNA expression levels of any genes in any tissues according to the GTEx portal (accessed on 30 January 2023) [53], suggesting that it is unlikely that these SNPs influence variations in opioid sensitivity among individuals by influencing the mRNA expression of some genes. The rs12652255 SNP was reported to be associated with the efficacy of Drotrecogin alfa, a drug with antithrombotic, profibrinolytic, anti-inflammatory, and cytoprotective properties in patients with severe sepsis [54], but its contribution to the efficacy of opioid analgesics remains unknown.
Among the candidate SNPs that were selected in our GWAS for the 0-12 h plasma MEC in the dominant model, the rs9533839 SNP was included, which was annotated as the SERP2 gene (Supplementary Table S2). This gene was also significantly associated with the same trait (Supplementary Table S4) and the 0-6 h plasma and 0-12 h effect-site MECs (Tables 3 and 4) in the gene-based analysis. The SERP2 gene encodes stressassociated endoplasmic reticulum protein family member 2 (SERP2), which is predicted to be involved in the endoplasmic reticulum unfolded protein response and protein glycosylation. Although SERP2 mRNA is known to be highly expressed in the brain, followed by the testis, according to the National Center for Biotechnology Information (NCBI) database, the functional relationship between this protein and the opioid system is unknown. In human cytogenetic studies, microdeletions of the SERP2 gene were reported to be associated with acute lymphoblastic leukemias in children with Down syndrome [55], and focal deletions of this gene were also identified in 2-6% of adult cases of acute lymphoblastic leukemia [56]. However, no genetic variants, such as SNPs, in this gene region have been reported to be associated with diseases or other phenotypes. HaploReg v. 4.1 showed that the rs9533839 SNP could change six motifs for DNA-binding proteins and overlaps with an enhancer in nine tissues, and this SNP was found to be significantly associated with mRNA expression levels of the TUSC8 gene in the prostate, breast (mammary tissue), and minor salivary gland according to the GTEx portal (accessed on 30 January 2023). The TUSC8 gene encodes a non-coding RNA (ncRNA), TUSC8, and this ncRNA reportedly enhances the cisplatin sensitivity of non-small-cell lung cancer cells by regulating vascular endothelial growth factor A (VEGFA) [57], although the involvement of this ncRNA in opioid sensitivity remains unknown.
In the gene-set analysis, several significant associations were also found ( Supplementary Tables S6-S9). Some of the genes that were included in the gene sets were included in two or three kinds of gene sets (Supplementary Table S10). Among the gene sets that were included in two kinds of gene sets, the VEGFA gene (Supplementary  Table S10) is notable because VEGF-A protein, which is encoded by the VEGFA gene, is known to be involved in angiogenesis through activation of the opioid system [58], although opioids could also exert a proangiogenic effect at low doses but an antiangiogenic (toxic) effect at high doses [59]. Intriguingly, the ANGPT1 gene was included in the "pid_shp2_pathway" gene set, which was significantly associated with the 0-12 h effect-site MEC in the recessive model (Supplementary Table S9). The ANGPT1 gene encodes angiopoietin-1, a secreted glycoprotein that is a member of the angiopoietin family. Angiopoietin-1 is also known to be involved in angiogenesis. Mice that were engineered to lack angiopoietin-1 exhibited angiogenic deficits [60]. Although more studies are required, angiogenesis, with the involvement of angiopoietin-1, could also be modulated by actions of opioids, and the rs1283671 and rs1283720 SNPs within this gene region were found to be significantly associated with average daily opioid requirements for the treatment of cancer pain in our previous GWAS [33]. Excluded were patients with severe coexisting disease (ASA-PS ≥ III), those taking pain medication for chronic pain, and those who were unlikely to be able to use a PCA pump (e.g., because of dementia). All of the individuals who were included in the study were of Japanese origin. Peripheral blood samples were collected from these subjects for gene analysis. Detailed demographic and clinical data of the subjects are provided in Supplementary Table S1 and our previous reports [22,25].

Materials and Methods
The study was conducted according to guidelines of the Declaration of Helsinki and approved by the Institutional Review Board or Ethics Committee of Saitama Medical University International Medical Center and Tokyo Metropolitan Institute of Medical Science (Tokyo, Japan). Written informed consent was obtained from all of the patients.

Surgical Protocol and Clinical Data
The protocols for anesthesia, surgery, and postoperative pain management and clinical data are detailed in our previous reports [22,25]. Briefly, general anesthesia was induced with fentanyl (0.1 mg), propofol (1-2 mg/kg), and rocuronium (0.8 mg/kg). After tracheal intubation, the inhalation of sevoflurane (1.5% in inspired concentration) and continuous infusion of remifentanil (0.25 µg/kg/min) were started. General anesthesia was thus maintained with sevoflurane, remifentanil, and rocuronium. At the end of surgery, remifentanil and sevoflurane were discontinued, and fentanyl (usually ≥ 0.1 mg) was given for immediate postoperative pain relief. The average remifentanil infusion rate (in µg/kg/min) during surgery was calculated by dividing the total dose of remifentanil that was required during surgery by the duration of surgery and body weight. When patients complained of even mild abdominal pain, fentanyl was given in increments of 0.05 mg until sufficient pain relief was achieved.
Postoperative pain was then managed with i.v. fentanyl PCA using a PCA pump (CADD-Legacy Model 6300, Smiths Medical Japan, Tokyo, Japan) that was filled with 1000 µg fentanyl diluted with normal saline to a total volume of 100 mL. The demand dose, dose lockout time, maximum allowable demand dose per hour, and continuous rate were set at 20 µg (2 mL), 5 min, 12 times (240 µg), and zero, respectively. Patientcontrolled analgesia was principally continued for 24 h postoperatively. In cases of inadequate analgesia, i.v. flurbiprofen axetil (50 mg) or pentazocine (30 mg) was administered as a rescue analgesic. Severe postoperative nausea and vomiting were treated with i.v. droperidol (2.5 mg) or metoclopramide (10 mg).
Postoperative pain at rest was assessed on an 11-point numerical rating scale (0, no pain; 10, the worst pain imaginable). Sedation was assessed on a 4-point scale (0, awake and alert; 1, drowsy; 2, mostly asleep but easy to rouse; 3, asleep and difficult to rouse). Postoperative nausea and vomiting were assessed on a 4-point scale (0, no nausea or vomiting; 1, mild nausea; 2, severe nausea; 3, retching or vomiting). Postoperative pain scores, sedation scores, postoperative nausea and vomiting scores, respiratory rates, the cumulative number of PCA doses that were actually given to the patient, and the cumulative number of PCA doses that were attempted were recorded on a data collection sheet 2, 4, 6, 12, 18, and 24 h after surgery.
PCA fentanyl consumptions over 6 h, 12 h, and 24 h periods were calculated as cumulative doses of fentanyl that were actually given to patients via the PCA pump during the first 6 h, 12 h, and 24 h postoperative periods, respectively. The 6 h, 12 h, and 24 h total postoperative fentanyl requirements were calculated as sums of the i.v. fentanyl dose that was given around the end of surgery and 6 h, 12 h, and 24 h PCA fentanyl consumptions, respectively. Patient-controlled analgesia fentanyl consumptions and total postoperative fentanyl requirements were normalized to body weight. The 6 h, 12 h, and 24 h numbers of locked out doses were the differences between the cumulative number of doses attempted and doses that were given at 6, 12, and 24 h after surgery, respectively.
Plasma and effect-site concentrations of fentanyl over the 24 h postoperative period were estimated in each patient using BeConSim Monitoring (http://www.masuinet.com; accessed on 1 January 2020; Supplementary Figure S1)-a pharmacokinetic simulation program that was developed by Masui (2010) [40] (Supplementary Figure S4) based on Shafer's three compartment model [61]-by inputting relevant clinical data, including age, sex, height, weight, the fentanyl dose that was given around the end of surgery, and subsequent PCA fentanyl consumption profiles over the 24 h postoperative period. A pair of plasma and effect-site MECs of fentanyl were indicated by plasma and effectsite fentanyl concentrations that were estimated immediately before each self-dosing of fentanyl. All pairs of plasma and effect-site MECs in each patient were averaged over the 6 h, 12 h, and 24 h postoperative periods to determine pairs of average plasma and effect-site MECs over these periods, which were expressed as the 0-6 h, 0-12 h, and 0-24 h plasma and effect-site MECs, respectively. Because many patients had completely consumed PCA fentanyl (1000 µg) by 24 h postoperatively and not by 12 h postoperatively, 0-6 h and 0-12 h plasma and effect-site MECs and not 0-24 h plasma and effect-site MECs were used for the main study endpoints. Detailed clinical data of the subjects are detailed in Supplementary Table S1. For patients who underwent LAC, 10 mL of venous blood was sampled during anesthesia for the later preparation of genomic DNA specimens. After total genomic DNA was extracted from whole-blood samples using standard procedures and the concentration was adjusted to 100 ng/µL, whole-genome genotyping was performed using the Infinium Assay II with an iScan system (Illumina, San Diego, CA, USA) ac-cording to the manufacturer's instructions. A total of 921,239 SNP markers survived the entire quality control filtration process and were used for the GWAS (detailed in our previous report) [22]. Two kinds of BeadChips were used for genotyping 256 and 95 samples, respectively: HumanOmniExpressExome-8 v. 1.0 (total markers: 951,117) and HumanOmniExpressExome-8 v. 1.1 (total markers: 958,178). Approximately 926,000 SNP markers were commonly included in all of the BeadChips. Quality control was properly performed the same way as in our previous report [22,25].
For phenotypes of the estimated 0-6 h plasma MEC and 0-12 h effect-site MEC, log QQ p-value plots as a result of the GWAS for the combined 351 samples were subsequently drawn to check the pattern of the generated p-value distribution, in which the observed p values against the values that were expected from the null hypothesis of a uniform distribution, calculated as −log10 (p value), were plotted for each model. All of the plots were mostly concordant with the expected line (y = x), especially over the range of 0 < −log10 (p value) < 4 for each model, indicating no apparent population stratification of the samples that were used in the study (Supplementary Figures S1 and S2).

Gene-Based and Gene-Set Analyses
Gene-based and gene-set approaches were adopted with Multi-marker Analysis of GenoMic Annotation (MAGMA) v. 1.06 [43], which is also available on the Functional Mapping and Annotation of Genome-Wide Association Studies (FUMA GWAS) v. 1.3.9 platform [44], to better understand genetic backgrounds and molecular mechanisms that underlie complex traits, such as opioid sensitivity in patients who underwent LAC. To examine the combined relationship between all genetic markers in the gene and the phenotype, genetic marker data were aggregated to the level of full genes in the gene-based analysis. Similarly, individual genes were compiled into groupings of genes with similar biological, functional, or other properties for the gene-set analysis. Gene-set analyses can thus shed light on the role that particular biological pathways or cellular processes may play in the genetic basis of a trait [43]. In these analyses, associations were explored for genes on autosomes 1-22 and the X chromosome, and the window of the genes to assign SNPs was set to 20 kb, thereby assigning SNPs within the 20 kb window of the gene (both sides) to that gene. For the reference panel, the 1000 Genome Phase3 EAS population was selected (http://ftp.1000genomes.ebi.ac.uk/vol1/ftp/release/20130502; accessed on 18 January 2023). In the gene-set analysis, gene sets were defined using the Molecular Signatures Database (MSigDB) v. 7.5.1 (https://www.gsea-msigdb.org/gsea/msigdb; accessed on 18 January 2023) [62]. A total of 10,678 gene sets (curated gene sets: 4761, Gene Ontology [GO] terms: 5917) from MsigDB were tested. In both analyses, Bonferroni correction for multiple testing was performed for all tested genes and gene sets. Adjusted values of p < 0.05 in the results were considered significant. The FUMA GWAS platform was also used for the visualization of QQ plots for the GWAS results and Manhattan plots for the gene-based analysis results.

Statistical Analysis
A three-stage GWAS was conducted for patients who underwent LAC to investigate the association between opioid sensitivity after surgery and the 921,239 SNPs that met the quality control criteria in a total of 351 subjects (117, 117, and 117 subjects for the first-, second-, and final-stage analyses, respectively) for whom postoperative clinical data were available, as described in our previous report [22]. As an index of opioid sensitivity after surgery, the estimated 0-6 h plasma MEC, 0-12 h plasma MEC, 0-6 h effect-site MEC, and 0-12 h effect-site MEC were used because these calculated values were expected to reflect the efficacy of fentanyl in each individual. Prior to the analyses, the quantitative values (ng/mL) were natural-log-transformed for approximation to the normal distribution according to the following formula: Value for analyses = Ln (1 + MEC value [ng/mL]). To explore the association between the SNPs and phenotypes, linear regression analyses were conducted in each stage of the analysis, in which the MEC value (ng/mL; log-transformed) and the genotype data for each SNP were incorporated as dependent and independent variables, respectively, with covariates that were found to be strongly associated with the dependent variable in a preliminary study. Male genotypes were not included in the analysis of X chromosome markers, whereas both male and female individuals were included in the association study for autosomal markers. Additive, dominant, and recessive genetic models for each minor allele were used for the analyses because of the previously insufficient knowledge about genetic factors that are associated with opioid sensitivity. The GWAS procedure is summarized in Supplementary Figure S5 for the 0-6 h plasma MEC and 0-12 h effect-site MEC, and the procedure was similar for the 0-12 h plasma MEC and 0-6 h effect-site MEC (details not shown). In the first-stage analysis of 117 subjects, the SNPs that showed statistical p < 0.05 were selected as candidate SNPs for the second-stage analysis among the 921,239 SNPs. For these SNPs, the second-stage analysis was conducted, and SNPs that showed p < 0.05 for the single analysis of this stage and combined analysis of the first and second stages were considered possible candidates. Similarly, the final-stage analysis was conducted by setting the threshold p values at 0.05, in which the SNPs that showed p < 0.05 for the single analysis of this stage and combined analysis of the first, second, and final stages were considered possible candidates. The potent SNPs were selected from these SNPs after LD-based SNP pruning to remove redundant SNPs due to strong LD (threshold r 2 = 0.8) with each other, as conducted in a previous report [63]. In the final stage, q values of the false discovery rate (FDR) were also calculated to correct for multiple testing for the SNPs that were selected after the second-stage analysis and LD-based SNP pruning, based on previous reports [64,65]. The SNPs that showed q < 0.05 in the analysis among the SNPs that were selected after the final stage were considered to be genome-wide significant. Hardy-Weinberg equilibrium was additionally tested using Exact Tests for genotypic distributions of SNPs that were significantly associated with the phenotype. To calculate q values, Stratified False Discovery Rate (SFDR) software (http://www.utstat.toronto.edu/sun/Software/SFDR/index.html; accessed on 18 January 2023) was used [64][65][66]. All of the statistical analyses for genetics were performed using gPLINK v. 2.050, PLINK v. 1.07 (http://zzz.bwh.harvard.edu/plink/index.shtml; accessed on 18 January 2023) [67], and Haploview v. 4.2 (https://www.broadinstitute. org/haploview/haploview; accessed on 18 January 2023) [68]. Additionally, correlation analysis, the Mann-Whitney test, and linear regression analysis were performed for statistical analyses of clinical variables using SPSS Statistics v. 25 software (IBM, Armonk, NY, USA). For the statistical analyses of clinical variables, the criterion for significance was set at p < 0.05.

Power Analysis
Statistical power analyses were preliminarily performed using G*Power v. 3.0.5 [69] as previously described [22,25]. Power analyses for the linear regression analyses revealed that the expected power (1 minus type II error probability) was 98.6% for Cohen's conventional "medium" effect size of 0.15 [70] when the type I error probability was set at 0.05 and sample sizes were 117, corresponding to the sample size of each stage analysis in the present study. However, for the same type I error probability and sample sizes of 117, the expected power decreased to 32.9% when Cohen's conventional "small" effect size was 0.02. Conversely, the estimated effect sizes were 0.0682 for the same type I error probability and sample sizes of 117 to achieve 80% power. Therefore, a single analysis in the present study was expected to detect true associations with the phenotype with 80% statistical power for effect sizes from large to moderately small, but not too small, although the exact effect size has been poorly understood in cases of SNPs that significantly contribute to opioid sensitivity.

Linkage Disequilibrium Analysis
The LD analysis was performed using Haploview v. 4.2 [68] for a total of 351 samples from patients who underwent LAC for the genomic position from~174,760,000 tõ 174,900,000 on chromosome 5 (GRCh37) that includes both the rs966775 SNP and DRD1 gene and its flanking region to identify relationships between SNPs in this region. The commonly used D and r 2 values were pairwise calculated using the genotype dataset for each SNP to estimate the strength of LD between SNPs. Linkage disequilibrium blocks were defined as in a previous study [71]. For the visualization of LD plots with information about genomic position and related gene transcripts, the LDmatrix tool was also used (https://ldlink.nci.nih.gov/?tab=ldmatrix; accessed on 22 January 2023).

Reference of Databases
Several databases and bioinformatic tools were referenced to more thoroughly examine the candidate SNP that may be related to human opioid analgesic sensitivity, including the NCBI database (http://www.ncbi.nlm.nih.gov; accessed on 19 January 2023), HaploReg v. 4.1 (https://pubs.broadinstitute.org/mammals/haploreg/haploreg. php; accessed on 19 January 2023) [51], SNPinfo Web Server (https://snpinfo.niehs.nih. gov; accessed on 19 January 2023) [52], and Genotype-Tissue Expression (GTEx) portal (https://gtexportal.org/home/; accessed on 19 January 2023) [53]. HaploReg is a tool for investigating non-coding genomic annotations at variations in haplotype blocks, such as potential regulatory SNPs at disease-associated sites [51]. The SNPinfo Web Server is a set of web-based SNP selection tools (freely available at https://snpinfo.niehs.nih.gov; accessed on 19 January 2023) where investigators can specify genes or linkage regions and select SNPs based on GWAS results, LD, and predicted functional characteristics of both coding and non-coding SNPs [52]. The GTEx project, an ongoing effort to create a comprehensive public resource to research tissue-specific gene expression and regulation [53], is the basis for the GTEx portal, which offers open access to such data as gene expression, quantitative trait loci, and histology images.

Conclusions
In conclusion, our GWASs revealed that the rs966775 SNP and SERP2 gene were significantly associated with estimated plasma MECs over the 0-6 h and 0-12 h postoperative periods of fentanyl that was administered for the treatment of postoperative pain in LAC patients. Although the present results need to be corroborated by more research with larger sample sizes and greater statistical power, these findings indicate that the rs966775 SNP near the DRD1 and SERP2 genes could serve as markers that predict the efficacy of opioid analgesics for the treatment of postoperative pain.  Informed Consent Statement: Informed consent was obtained from all subjects involved in the study.
Data Availability Statement: Data that are presented in this study are available upon request from the corresponding author.